Phase transition between synchronous and asynchronous updating algorithms 
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We update a one-dimensional chain of Ising spins of length L with algorithms which are param- 
eterized by the probabihty p for a certain site to get updated in one time step. The result of the 
update event itself is determined by the energy change due to the local change in the configura- 
tion. In this way we interpolate between the Metropolis algorithm at zero temperature for p of the 
order of 1/L and for large L, and a synchronous deterministic updating procedure for p = 1. As 
function of p we observe a phase transition between the stationary states to which the algorithm 
drives the system. These are non-absorbing stationary states with antiferromagnetic domains for 
p > Pc and absorbing states with ferromagnetic domains for p < Pc- This means that above this 
transition the stationary states have lost any remnants to the ferromagnetic Ising interaction. A 
measurement of the critical exponents shows that this transition belongs to the universality class of 
parity conservation. 

PACS numbers: 05.70.Ln , 05.50.+q , 64.90.-|-b 
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The issue of synchronous versus asynchronous up- 
dating algorithms has attracted much attention in 
connection with Boolean networks [E S S S) IB> 
ral networks biological networks and game 

theory [lo| . In a synchronous updating scheme all 
the units of the system are updated at the same time. 
Asynchronous updating usually means that not all units 
are updated at the same time. The algorithm can be 
asynchronous in the sense that each unit is updated 
according to its own clock, as in distributed systems 
for parallel processing [ll[, or it is asynchronous in the 
sense that only one randomly chosen unit is updated 
at each step, as in Monte Carlo algorithms [T5|. In the 
context of thermal equilibrium dynamics updating algo- 
rithms like the Metropolis algorithm are designed 
in a way that they drive the configurations to a set 
that is representative for the Boltzmann equilibrium 
distribution. The dynamics of the algorithm then enters 
only in an intermediate step, it is not representative for 
the intrinsic dynamics of the system that is determined 
by the Hamiltonian. In out-of-cquilibrium systems the 
updating scheme may play a more prominent role. The 
number of attractors in Boolean networks, for example, 
increases exponentially with the system size IJ] for syn- 
chronous update, and with a power for critical Boolean 
networks [15] for asynchronous update. The phase 
diagrams of the Hopfield neural network model [y, [3| 
and the Blume-Emery-Griffiths model [13, depend 
on the updating mode as well, while those of the Q-state 
Ising model [19| and the Sherrington-Kirkpatrick spin 
glass [20'| are independent on the used scheme. 
Probably neither a completely synchronous nor a random 
asynchronous update is realistic for natural systems. 
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Here we interpolate between these two extreme cases 
not in a more realistic way, but in a way that allows to 
identify a phase transition between the stationary states. 
We focus our attention on a one-dimensional Ising model 
at zero temperature. We visit all the sites and select each 
of them with probability p. Then we update simultane- 
ously each of the selected candidates by flipping their 
spins if the local energy of the system is not increased 
due to this change. Tuning the value of p, we are able 
to interpolate the algorithm from an asynchronous one, 
in the thermodynamical limit corresponding to the 
Metropolis algorithm, to a synchronous one. We observe 
a phase transition between the stationary states from 
non-absorbing for p > pc to absorbing ones for p < Pc, 
with Pc the critical threshold of the transition. The 
time evolution of the Ising chain can be represented as a 
directed percolation (DP) process if bonds between spins 
of opposite signs are called active. As our measurements 
of the critical exponents show, the transition between 
these stationary states belongs to the universality class 
of parity conservation (PC). 

We consider a one-dimensional lattice of length L. 
To each site i, i — 1,...,L of this chain we assign 
a spin variable ct^, where ai takes the values -1-1 
or —1. The Hamiltonian of the system is given by 
H = — J ^f^-^^ cTiai+i , where J is the coupling con- 
stant between neighboring sites. Here we consider 
ferromagnetic couplings so that J = 1. Moreover, for 
definiteness we choose periodic boundary conditions so 
that (Tl+1 = o'l- The results will not depend on this 
choice. 

A standard local stochastic updating algorithm for 
studying thermodynamic equilibrium of the Ising model 
is the Metropolis algorithm l^, [13) IHj. Given a 
configuration of the system I](t) = {o'i{t)} at time 
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t, we pick up one site j at random and flip its spin 
with probability Pj{t) = min {1, exp [— /JAE'j (i)]}. Here 
f3 = 1/kT, with k the Boltzmann constant and T the 
temperature, while AEj{t) = 2(Tj{t) [aj-i{t) + <7j^i{t)] 
is the difference in energy that a flip of aj would induce. 
In particular, for zero temperature an increase in energy 
of the resulting configuration is always rejected. After 
the single update of the j-th site the time increases by 
t^t+l/L~t' (and L single updates are considered 
as one time unit). The new configuration is given by 
S (i') = {ai{t);aj{t')}, where all sites i ^ j have the 
same spin value as they have at time t, but only the 
spin of the j-th site is eventually flipped. It is well 
known that the Metropolis algorithm drives our chain 
of Ising spins to the Boltzmann equilibrium distribution 
for a sufficiently large number of update events. In 
the context here we emphasize that the Metropolis 
algorithm is fully asynchronous in the sense that we 
have only one spin flip per single update event. 
Here we are no longer interested in the equilibrium 
properties of the Ising model, but in the algorith- 
mic dynamics applied to Ising spins. Therefore we 
give up the asynchrony of the algorithm. Given the 
conflguration E(t) = Wi{t)} Sit time t, we visit all 
sites and select each of them with probability p. The 
selected sites are ji, . . . , jm and m, the total number of 
selected sites, is a random integer obeying the binomial 
distribution B{m,L,p) = (^)p"' (1 - p)^"". Each 
of the m selected sites is then updated according to 
the Metropolis rule at temperature zero, so that the 
spin of the j„-th site is flipped with probability Pj^^t), 
V u = 1, . . . , m. After one step of the algorithm, the time 
increases as t ^ t + p = t' , and the new configuration 
is S(t') = {ai{t);aj^{t');...;aj^{t')}, where all sites 
i ^ ji, ■ ■ ■ ,jm, have the same spin value as they have 
at time t, while the spins of the m selected sites are 
eventually fiipped. One time unit has passed when the 
average number of update events equals to the total 
number of sites L. By varying p we can "tune" the 
algorithm from asynchrony, for p of order of 1/L, to 
synchrony, for p = 1. In particular, for p of the order 
oi 1/L and for sufficiently large values of L we recover 
the usual Metropolis algorithm when two or more 
simultaneous selections which occur with probability 
at most of the order of become negligible. We 

stress the fact that the spin value of the selected sites 
Ji ) • ■ ■ J jm at time t' is eventually flipped according to 
the actual values of the spins at time t, so that they 
depend only on the configuration Differently from 

the standard Metropolis algorithm at zero temperature, 
for general values of p, the total energy of the new 
configuration S(t') can be increased with respect to the 
total energy of the old configuration Moreover, the 

long-time configurations — > oo) do no longer obey 
the Boltzmann distribution. 

In this paper, for simplicity, we focus on the case of 



zero temperature. The equilibrium ground-state of the 
one-dimensional Ising ferromagnet at zero temperature 
is one of the two ferromagnetic states with all spins 
positive or negative. In contrast, the completely syn- 
chronized dynamics does not drive the system to the 
ground state, but acts as parallel algorithm and amounts 
to a deterministic map T: E(t + 1) = TE(i), for all t. 
After a transient time Iq < L/2, the algorithm drives 
the system into a cycle of length two [19|, where the 
system "flips" between two configurations, Ei = {^(Ji} 
and E2 — {^(Ti}, such that E2 = TEi and Ei = TE2, 
for all t > to. In particular, these conflgurations result 
from each other by an overall flip of signs in the sense 
that E2 = {— ^cTi} and Ei = { — ^CTi}. Moreover, 
these states are anti-ferromagnetic because we observe 
domains (neighboring sites with the same value of the 
spin) of a length of at most two sites. Therefore, it is 
natural to study intermediate values of p, in particular 
to focus on the transition between the ferromagnetic and 
the anti-ferromagnetic conflgurations of the flnal state. 
Let us consider the active bonds of the system, where 
we define a bond as active if it connects two sites with 
opposite spins. As remnant of the zero-temperature Ising 
model only sites belonging to at least one active bond 
can flip and do flip if they are selected as candidates 
for the updating. Only a few elementary processes, 
that involve active bonds, can take place: diffusion, 
annihilation and creation. For clarity of notation, let 
us indicate as t a site with positive spin and as | a 
site with negative spin. Consider, for example, a local 
configuration such as • • • ttii • • • at time t: at time t+p 
it can evolve to • • • ttt| • • • or to ■ • • ^iU • • • , depending 
on whether they are selected for update that happens 
with probability 2p{l — p) [diffusion], or to • • • titi ■ ■ • 
with probability p^ [creation], or it remains unchanged 
with probability (1 —p)^- Using the same rules, a local 
configuration such as • • • TTiTT ■ ' ■ at time t later, at 
time t+p, can become • • • titit ' ' ' with probability 
p^ [creation], or ■•• ... with probability p{l — p)"^ 

[annihilation], etc. ... It is easily checked that the parity 
of active bonds, that is the number of active bonds 
modulo 2, is conserved. 

The former considerations suggest that the transition 
between ferromagnetic and anti-ferromagnetic behavior 
(without active bonds and with dominance of active 
bonds, respectively) is given by the competition of 
annihilation and creation of active bonds. In particular, 
the creation is favored by a synchronous updating 
scheme, because a new couple of active bonds can be 
created only if two neighboring sites simultaneously flip 
their spins. Therefore the transition between stationary 
states (ferromagnetic and anti-ferromagnetic ones) 
can be considered as a transition between the asyn- 
chronous/synchronous updating schemes. Furthermore, 
this transition corresponds to a (1 + l)-dimensional 
DP transition 23). The active bonds correspond 
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cal, critical and supercritical regimes. We show only the part ^q2 ^^a " .|g3 

around the initial active bond. From left to right: p — 0.39, |_ |_ 

p = 0.41 =p^,p = 0.43. 



to the occupied sites in DP. A qualitative picture of 
this transition is shown in Figure [TJ We plot the time 
evolution of an isolated active bond (that is not the 
only one in the system) for three different values of p. 
We only display the intermediate part of the lattice 
around the initial active bond. The initial configuration 
is chosen as ■ • ■ ttii ■ ' ■ • In the subcritical regime 
p < Pc, the initial active bond will remain alone for 
most of the time, diffusing, and, from time to time, 
creating couples of new active bonds which annihilate 
soon. In the supercritical regime p > pc, the average 
number of creations is larger than the average number of 
annihilations, so that the active bonds spread over the 
entire system. In the critical regime p = Pc, annihilation 
and creation processes are balanced and the active bonds 
do not spread over the system, but remain confined in a 
finite region leading to the same qualitative picture as 
we see in the subcritical region. 

In order to obtain a quantitative description of the 
transition, we use as order parameter the density of 
active bonds p, given by the ratio of the number of active 
bonds and the total number of bonds in the lattice. The 
initial condition is always chosen as p(0) = 1, so that 
the lattice is fully occupied with active bonds. All data 
points are obtained from averaging over at least 10'^ 
realizations and up to 10^ realizations for small sizes 
of the lattice. Here, these values of L are considered 
as the large-volume limit in the following. Let us first 
determine the critical probability pc- In Figure HJ'V) we 
plot the time behavior of p for three different values of 
p and L = 10^. As we can see, for p — 0.41 p{t) ~ 
with 6 = 0.286(1), while for p — 0.40 p{t) decreases with 
negative curvature, for p = 0.42 p(t) it increases with 
positive curvature. The positive curvature characterizes 
the different phase. Therefore, within the given accuracy, 
we locate the critical threshold Pc as the largest value 
such that the curvature is non-positive. In this way we 
obtain pc — 0.41(1) as the critical point. Determining 
Pc with higher precision would amount to an increase in 
the linear size of the lattice L and the time for observing 



Figure 2: A) Time decay of the density of active bonds p 
in the subcritical {p = 0.40), critical (p — 0.41 — pc) and 
supercritical {p = 0.42) regime, from bottom to top, respec- 
tively. B) Average density of active bonds ps surviving up 
to time t over samples, plotted as a function of time at the 
critical point for L = 16, 32, 64, 128, 256, 512, 1024, 2048, 4096 
and 6144, from top to bottom respectively. The dotted line 
has a slope of —S = —0.286(1). C) Average density of active 
sites over the surviving samples pa as a function of size L of 
the lattice for p = 0.40 (empty dots) 0.41 (fuU dots) and 0.42 
(empty squares). The full line has slope —j3/v±^ = —0.51(1). 
D) Average relaxation time t as a function of the size L of 
the lattice for p = 0.40, 0.41 and 0.42 [the symbols are the 
same as in C)]. The full line has slope z = i'\\/i'± = 1.746(2). 



the positive curvature of p. Here we do not perform 
this kind of computationally expensive simulations, but 
calculate the critical exponents without increasing the 
precision in pc- 

Critical exponents. Along the determination of Pc 
we have already read off the exponent i5 from the time 
evolution of p: for p < pc and in the thermodynamic 
limit L oo, we should observe p{t) ~ t~°'^P\ with 
a{p) a continuous and monotonic decreasing function 
of p and enclosed by curves with a{p) = 1/2 for 
p 0+, as in the case of the standard Metropolis 
algorithm [iT 



2l[, and a(j>c) — S, as observed in our 



numerical simulation. Therefore, the negative curvature 
observed in the subcritical regime, is a finite-size effect. 
Also in the supercritical regime the curve bends down to 
zero as a finite-size effect after a sufficiently large time, 
but the positive curvature signals the onset of the new 
(supercritical) phase. 

Finite-size scaling analysis. Finite-size effects in the 
density of active bonds are manifest in two ways: if we 
follow the time evolution of a certain configuration of 
a chain of length L, we observe a power-law decay of 
active bonds up to a certain time td{L,p) for p < pc- 
After that, either the density drops to zero faster than 
a power-law, this happens for most configurations, for 
which a fluctuation drives the system into the absorbing 
state, or, in the minority of evolutions, the number of 
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active bonds fluctuates around a plateau, before the 
plateau drops to zero in the end. Configurations of this 
minority are called surviving configurations up to time 
T. Now it is easier to locate the onset of the plateau 
than the onset of a faster decay, and therefore to study 
the finite-size scaling of the density of active bonds ps 
of surviving configurations as a function of L, averaged 
only over the surviving realizations 0| . In Figure [2j3) 
we plot the time behavior of ps at the critical point for 
several values of the size L. As we can see, after an 
initial transient in which ps decreases as , it reaches a 
stationary value depending on L. This value vanishes in 
the thermodynamical limit, since the plateau is a finite- 
size effect. We average along the values of the plateau 
(since averaging over time and over different realizations 
are assumed to be equivalent), for different system sizes. 
In Figure [2p) this average value is plotted as a function 
of L for three values of p {p — 0.40, p — 0.41 = pc 
and p = 0.42). Again, at the critical point wc find a 
power law decay ps L~^/'^^, with f3/iy±_ ~ 0.51(1), 
while the decay deviates from the power-law behavior in 
both the subcritical and the supercritical regimes. (Here 
the exponent (3 characterizes the behavior of the order 
parameter p{p) ^ \pc — p\^ , the spatial correlation 
length £_{p) ^ \pc — p\~'^^ both in the vicinity of Pc-) 
Moreover, from the finite-size scaling analysis we 
calculate the dynamical exponent z — vu/i/^, where 



'II 



1/ 

characterizes the time-like correlation length 
t(p) ~ \pc — p\~'^K The exponent z is derived from the 
relaxation time r needed by a finite system to reach the 
absorbing configuration [p{t > r) = 0]. In Figure [2jD) r 
is plotted as a function of L ior p = 0.40, p = 0.41 = pc 
and p — 0.42. At p — pc we find again a power-law 
dependence r ~ L^, with z — v\\/v± — 1.746(2), while 
for different values of p r behaves differently from a 
power- law behavior, as it is seen in Figure [2jl)). Further- 
more, in Figure[3]we verify the finite-size scaling relation 
p{pc,L,t) - L-I^^^^f (t/L'^), where /(•) is a suitable 
universal function. Finally, we determine the static 
exponent 7 from the growth of fluctuations in the order- 
parameter susceptibility, defined as x = L{{p'^) — (p)^), 
via Xs ^ L^/'^^ at (where the index s again refers 
to an average over the surviving configurations). We 
find 7 ~ 0, but we do not report any figure here. 
We further determine the dynamical exponent 77 and 
check the exponents 5 and z in the following way. 
Starting from a configuration like ••• ttiTT '''i that 
is a ferromagnetic configuration with only two active 
bonds, we numerically compute the survival probability 
P{t) (that is the probability that the system had not 
entered the absorbing state up to time t), the average 
number of active bonds n{t) and the average mean 
square distance of spreading R?{t) from an arbitrary 
selected site. These quantities are expected to behave 
at the critical point Pr as P{t) ~ , n{t) ^ and 
i?2(t) ^ i2/z [22,[2|,[2|]. From these simulations we find 




t/L^ 



Figure 3: The main plot shows the finite-size scaling of the 
density of active bonds p. The inset shows the unsealed data, 
where the dotted fine has a slope equal to —5. The datasets 
are the same as in Figure [2] 



77 2± and values of 6 and z consistent with the former 
ones within the error bars. The figures are not displayed. 

In conclusion, we have calculated the static expo- 
nents P/i'± = 0.51(1) and 7 ~ and the dynamical 
exponents S = 0.286(1), = z = 1.746(2) 

and 77 ~ 0. These values arc consistent with 
the ones conjectured for the PC universality class 
f3/iyj_ = 1/2,7 = 0, (5 = 2£7,iA\/i^ = 21/12 and 

77 = [22, lii, [ii [ia S [B, 0, S [3^. This is 

expected because periodic boundary conditions for Ising 
spins as well as the updating rules we used preserve 
the parity of the number of active bonds, and for free 
boundary conditions only the boundaries may violate 
the conservation of parity which leads to a negligible 
effect. Even if parity is non-conserved as in case of the 
non-equilibrium kinetic Ising model (29l . [30| , for which 
spin-fiip and spin-exchange dynamics are mixed, the 
PC universality class is observed. An increase in the 
dimensionality of the lattice would allow a mean-field 
description of the transition, since the PC universality 
class has the critical dimension dn — 2 22. \2c 



Oi 



preliminary numerical simulations for this case show 
qualitatively similar phases as in d = 1, but with the 
transition point shifted towards 1. Our work represents 
the first systematic interpolation between a synchronous 
and asynchronous updating scheme. The existence of a 
phase transition sheds some light on the interpretation of 
stationary states whenever they depend on the updating 
mode. Beyond the transition these states may have 
lost any remnants of the "intrinsic" dynamics (in our 
case the ferromagnetic Ising interaction). Instead, they 
are representative for the dynamics of the updating 
algorithm itself. 
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